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Abstract. Large-scale reionization simulations are described which combine the results of cos- 
mological N-body simulations that model the evolving density and velocity fields and identify the 
galactic halo sources, with ray-tracing radiative transfer calculations which model the nonequilib- 
rium ionization of the intergalactic medium. These simulations have been used to predict some of 
the signature effects of reionization on cosmic radiation backgrounds, including the CMB, near-IR, 
and redshifted 21cm backgrounds. We summarize some of our recent progress in this work, and 
address the question of whether observations of such signature effects can be used to distinguish the 
relative contributions of galaxies of different masses to reionization. 
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1. INTRODUCTION 



The talk on which this paper is based summarized some of our recent progress towards 
simulating cosmic reionization and the radiation backgrounds from this epoch 2 . Here we 
are constrained by length limitations to be even more selective, so we will not presume 
to review this large and active field, but rather will focus on just a few of our most recent 
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results mentioned in the talk. 

When the first stars and galaxies formed, they released ionizing radiation into the 
surrounding intergalactic medium (IGM), creating a patchwork quilt of H II regions that 
grew in size and number and merged over time, until their overlap eventually ionized the 
entire IGM. The spectra of distant quasars and galaxies seen through this IGM suggest 
this overlap occurred by redshift z ov > 6, while fluctuations in the polarization of the 
CMB on large angular scales imply that reionization must have begun at significantly 
higher redshift, to provide enough intergalactic, electron-scattering optical depth, x es , 
to explain this. Observing cosmic reionization, then, would seem to be a good way to 
learn about first star and galaxy formation. While direct observation of the sources of 
reionization remains a challenge, their impact on cosmic background radiation during 
the epoch of reionization (EOR) offers another probe. 

Predictions require the modelling of the inhomogeneity of reionization on very large 
scales (> 100/h cMpc), starting from realistic cosmological initial conditions, in order 
to be statistically meaningful. Ideally, such models would simultaneously resolve the 
formation and internal structure of individual galaxies and, inside each, the births and 
deaths of all the stars, including exchange of mass, energy and radiation between the 
IGM and the galaxies and the back-reaction caused by reionization, itself. The dynamic 
range required for this "ideal" is currently very far from achievable, however. 

Fortunately, there is some separation of scales which is possible as a first approxi- 
mation, in which the location sites of reionization sources are identified with the dark- 
matter-dominated galactic halos that host galaxies, their efficiency as sources is param- 
eterized in some phenomenological way, and their radiation is transferred through the 
evolving IGM, whose changing ionization state also changes its opacity. This approxi- 
mation is aided by the fact that reionization proceeds by propagating ionization fronts 
(I-fronts) that travel highly supersonically through the IGM, the I-front regime known as 
"weak, R-type," so the H II region boundaries race ahead of the hydrodynamical back- 
reaction of the gas in response to its pressure increase by photoheating [1]. While this 
back-reaction can affect small-scale structure, it does not directly affect the large-scale 
pattern of reionization patchiness. On the other hand, it does affect the efficiency of 
small-mass galaxies as sources, since the pressure increase of the photoheated IGM op- 
poses baryonic gravitational collapse on small-scales, sometimes referred to as "Jeans- 
mass filtering" [2, 3], so small-mass halos inside IGM H II regions may have signif- 
icantly reduced star formation rates. By changing the source efficiencies in this way, 
such feedback can alter the large-scale structure of reionization, too. The mass upper- 
limit of halo sources suppressed this way is uncertain, but estimates range at least to 
1O 9 M . The smallest galactic halos that can form stars, minihalos with mass between 
10 5 and 10 8 M (with virial temperatures T < 10 4 K), are especially vulnerable to such 
radiative suppression inside H II regions, since they can also lose the interstellar gas they 
already have, if overtaken by an 1-front as it crosses the IGM, which expels their gas by 
photoevaporation [4]. Minihalo star formation can also be suppressed by the rising UV 
background below 13.6 eV emitted by the same reionization sources, of H2-dissociating 
radiation in the Lyman-Werner bands (e.g. [5]), about which we will say more below, 
which limits their contribution to reionization even in neutral patches of the IGM. Mini- 
halos depend upon H2 molecules to cool their gas below their virial temperatures to 
enable star formation, since atomic cooling in primordial composition gas is ineffective 



below 10 4 K. 

In what follows, therefore, we will group the sources of reionization according to 
their halo mass, in three broad bins: MiniHalos (MH), Low-Mass Atomic-Cooling Halos 
(LMACH) with 10 8 <M < 10 9 M , and High-Mass Atomic-Cooling Halos (HMACH) 
above 1O 9 M . MH and LMACH sources are both suppressed inside H II regions, 
MHs are suppressed in neutral zones, too, if the mean intensity of LW background 
radiation there exceeds some threshold value, 7lw > ^LW,th> an d HMACHs are immune 
to suppression. We will say more about the relative contributions of these three source 
types to reionization below, but first we describe the simulations we will apply to address 
that question. 

The numerical methods we have developed for this purpose are described in Section 
2. If this were a review of the field of reionization simulation, I would also mention the 
large body of work by other authors, too, that has gone into developing cosmological 
radiative transfer methods and applying them to this problem, on scales large and small, 
apart from our own, including [6-35]. To test these methods, we organized the Cos- 
mological Radiative Transfer Comparison Project several years ago, with workshops 
leading to published test problems and comparisons which have encouraged this rapid 
growth [36, 37]. 

2. REIONIZATION SIMULATIONS: N-BODY + RADIATIVE 

TRANSFER 

2.1. N-body Simulations 

For our second generation of reionization simulations, the density and velocity fields 
of the IGM and the source halos are given by large-scale cosmological N-body simula- 
tions by the new N-body code CUBEP 3 M, a P 3 M method, massively paralleled (MPI + 
Open MP) ([38, 39], J. Harnois-Deraps et al. 2012, in preparation). We mention three 
new simulations here with large particle number N, on a grid of N ce u s , in comoving pe- 
riodic box of size Lb ox > particle mass m p , and minimum resolved halo mass M m i n , given 
by (N, 7V cell , L box , m p , M min ) = (3073 3 , 6, 144 3 , 163 Mpc, 5 x 1O 6 M , 1O 8 M ), (5488 3 , 
10,976 3 , 29 Mpc, 5 x 1O 3 M , 1O 5 M ), and (5488 3 , 10,976 3 , 607 Mpc, 5 x 1O 7 M , 
1O 9 M ), respectively. Halos are identified "on-the-fly" by a spherical overdensity halo 
finder. Of these, the 163 (= 1 14/h) Mpc simulation with 29 billion particles resolves all 
LMACHS and HMACHs (i.e. > 1O 8 M ), of which there are more than 10 7 by z ~ 8. 
The 29 (= 20/h) Mpc simulation resolves MHs (> 1O 5 M ), the first of which form at 
z = 43, with more than 10 8 by z = 8, with length resolution as small as 182 pc. The 607 
(= 425 /h) Mpc box resolves all the HMACHS (> 1O 9 M ), the first forming at z = 26, 
with 4 x 10 7 halos by z = 8, in a volume comparable to that of the LOFAR EOR 21cm 
background survey. 



2.2. Radiative Transfer 



Radiative transfer (RT) simulations evolve the radiation field and nonequilibrium 
ionization state of the gas, in the density field provided by the N-body results above, 
smoothed to a coarser grid of ~ 256 3 to 512 3 cells, for different resolution runs, by 
ray-tracing the ionizing radiation from every galaxy halo in the box. A new, fast and 
efficient method was developed, the C 2 -Ray code (Conservative, Causal Ray-Tracing), 
which uses short characteristics and has the advantage that it accurately tracks I-fronts 
even using relatively large cell sizes (i.e. very optically thick when neutral) and time 
steps (i.e. much larger than the light-crossing and I-front crossing times per cell). While 
our first generation of reionization simulations used this method coded in Open MP 
for shared memory parallel computers [38, 40, 41], the second generation simulations 
highlighted here used a new version of C 2 -Ray, recoded for massively parallel computers 
with distributed memory (MPI + Open MP) [42-44]. The original C 2 -Ray only included 
hydrogen, but we recently generalized it to take account of helium and hard- spectrum 
sources (including X-rays), as well [45], and applied this to study the effect of adding 
quasar sources to a simulation which is otherwise reionized by starlight [46]. To track 
the inhomogeneous, H2-dissociating, LW background (1 1 .2 - 13.6 eV) from reionization 
sources, important for its feedback effect on MHs, we developed another algorithm 
which sums the contributions at each cell from all sources along the past light cone 
of that cell in every direction, taking account of cosmological redshifting, geometric 
dilution, and attenuation by the intervening H Lyman series opacity of the IGM [47, 48]. 

Every galaxy in the simulation emits ionizing radiation (unless suppressed) with ef- 
ficiency parameterized as follows. Within each halo mass group, we assume a constant 
M/L for simplicity, with fy ionizing photons released by each galaxy per halo baryon, 
where fy = f*f esc Ni, and /* is the star-forming fraction of halo baryons, / e sc is the ioniz- 
ing photon escape fraction, and Ni is the number of ionizing photons emitted per stellar 
baryon over the stellar lifetime. For example, if Ni = 50, 000 (top-heavy IMF), /* = 0.2, 
and /esc = 0.2, then f y = 2000, while if N t = 4, 000 (Salpeter IMF), /* = 0.1, and / esc = 
0.1, then fy = 40. This yields a source luminosity dNy/dt = /yMb ar y/(/im#A?*), for 
source lifetime (e.g. duration of burst) At* (e.g. 2xl0 7 years), M^y = halo baryonic mass 
= ^haio(^bary/^m), leading to a star formation rate SFR = (f Y / At*)(M\ }ai y/ f esc Ni) ~ 
(1.7M /yr)(/ y /4O)(O.l// esc )(4OOO/M)(lOMyr/A^)(M halo /lO 9 M ) [e.g. if f y = 40, 
/ esc = 0.1, /* = 0.1, and At* = 2;d0 7 , then SFR ~ (O.8M /yr)(M halo /lO 9 M )]. Since 
photons are released over the time interval At*, we define a second efficiency parameter 
which reflects this dependence, too: g y = fy/ (At*/ lOMyr) . 

When MH sources are included, we assume one massive Pop III star forms per halo 
(or an equivalent group of stars with the same effective /V,-value) before disrupting its 
own ISM and preventing further star formation, at least temporarily, in that MH, and 
/esc = 1- For such stars, we take At* = 1.92 Myr and Ni = 55,000, which implies an 
efficiency for halo mass Mh a i and stellar mass M*jjj given by fy : MH(M hSL \ , M*ju) 
= 338(M* ;/// /lOOM )/(M halo /lO 5 M ), where the average value of /^mh integrated 
over the halo mass function is comparable to the value at the low mass end of star- 
forming MHs. 

To account for the photoionization suppression of MHs and LMACHs, we set their 



fy to zero once the halos are inside H II regions. MHs are assumed suppressed by H2 
dissociation even in neutral zones, if they are located where 7 LW > J]_w,th, with value 
chosen in the range [0.01 — 0.1] x 10 _21 ergs _1 cm _2 sr _1 , found by other work as the 
threshold for suppressing MH star formation. 

In what follows, the reionization simulation cases are labelled to indicate the value of 
these parameters for different mass groups. The full notation that appears on some plots 
reads (L b ox)_g(^ r ,HMACH)_( ( ?yXMAC/f)S_M(M ;(!i ///)_J(7LW,th/10 _21 ), where "S" refers to 
LMACH suppression (e.g. 163Mpc_g8.7_130S_M300_J0.05). 

3. RESULTS 

3.1. Self-Regulated Reionization 

[41] demonstrated that, when LMACHs are included but suppressed inside H II re- 
gions, reionization is "self -regulated." The (more abundant) LMACHs start reionization 
earlier than the (less abundant) HMACHs, but the LMACH contribution saturates when 
only a fraction of the volume is ionized (i.e. long before z ov ), so reionization is finished 
by the HMACHs, which ultimately dominate reionization. This helped to explain how 
reionization might end late but extend in time enough to boost T es . The N-body simula- 
tions for this had to resolve all halos with M > 1O 8 M , and the C 2 -Ray simulations had 
to transfer the radiation from all of them, which limited the box size in practice to 50 
(= 35 /h) Mpc. To do this for a box large enough to make statistically meaningful pre- 
dictions of observables required us to develop the second generation of codes described 
in Section 2. Some of the results of these self -regulated reionization simulations for box 
size 163 Mpc and RT grid with 256 3 cells are described in [40, 43-45, 49-52]. 

To resolve all the individual MH sources, too, in so large a volume would require an 
N-body simulation with more than 30 trillion particles, however, well beyond current 
capabilities. Instead, we assign the MH contribution to each cell of the radiative transfer 
grid by a subgrid model based upon smaller-box N-body simulations like the 29 Mpc 
box mentioned above and a 9 Mpc box, both of which fully resolve the MH halo 
mass range, from which we derive the dependence of the MH mass function on local 
matter density, smoothed on the scale of the radiative transfer grid cells in our large-box 
reionization simulation [48] . With this, we have simulated reionization with MH sources, 
too, including their LW suppression, and found that a new kind of self-regulation occurs, 
in which MH sources start reionization even earlier than LMACHs, but they, too, saturate 
after ionizing only a fraction of the IGM. The later rise of the ACHs ultimately stops 
MH star formation altogether, and reionization is once again finished (and dominated) 
by the HMACHs, just as without MHs [48]. With MHs, however, reionization is greatly 
extended, which boosts T es and the large-angle polarization fluctuations of the CMB. 
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FIGURE 1. 21cm brightness temperature fluctuations from the EOR for Case L2 (see Table 1) 
(HMACH+LMACH, no MHs) for 2 box sizes (163 Mpc with 256 3 RT cells, 607 Mpc with 504 3 ). (Left 
panel) Rms fluctuations and skewness vs. frequency (dotted curve = 163 Mpc). (Right panel) Power 
spectrum vs. wavenumber at z > 7, when mean ionized fraction is x\ = 0.77 (leftmost curve = 607 Mpc). 



3.2. Very Large-Scale Reionization Simulations: 607 Mpc box 

To simulate self-regulated reionization in an even larger volume, to test if our simu- 
lation results for the 163 Mpc box have converged yet as a function of increasing box 
size and enable us to predict the 21cm background for a volume comparable to that of 
the LOFAR EOR survey, we have simulated reionization in the 607 Mpc box described 
in Section 2, with 504 3 RT cells. Our N-body simulation directly resolves all HMACHs 
(i.e. > 10 9 M Q ), but we also include the LMACHs by a similar subgrid prescription to 
that mentioned above for the MHs, by using smaller-box simulations with enough mass 
resolution to resolve LMACHs directly, to derive the LMACH mass function in each 
radiative transfer grid cell of the 607 Mpc box, as it varies with the cell-averaged matter 
density there. Some results have already appeared in [46], in which 21cm brightness 
temperature fluctuation maps were made for two cases, one ionized only by starlight, 
while the other included a quasar source, as well, to see if a matched-filter technique 
applied to the 21cm observations could successfully measure the presence and size of 
the H II region surrounding a quasar during the EOR. 

One comparison of interest of the 607 and 163 Mpc box simulations of self-regulated 
reionization is the 21cm brightness temperature predictions for each. In the high 
spin temperature limit generally assumed for this EOR signal, the differential bright- 
ness temperature 8T h = T b - T CM b = (28.74mK> H i(l + 5)[(1 +z)/10] 1 / 2 {l + [(1 + 
z) / H (z)](d\\\/ dr\\)}~ 1 , where xui is the neutral fraction, 8 is the overdensity, and vy is 
peculiar velocity component along the line of sight, observed at the redshifted frequency 
which takes peculiar velocity into account. In Figure 1 (left panel), we compare the 



TABLE 1. Reionization simulation results for the kSZ effect: global reionization history and 
CMB fluctuations. 



Case label* 


zov 


z 99% ~ z 20% 


z 75%~ z 25% 


, post-reion^ 


nreion 

^=3000 


ntotal 

^=3000 


LI 


8.3 


2.9 


1.8 


1.94 iiK 2 


0.83 juK 2 


2.77 /xK 2 


L2 


6.7 


1.8 


1.2 


1.69 juK 2 


0.66 jiiK 2 


2.35 jiiK 2 


L2M1 


6.7 


6.6 


1.8 


1.69 juK 2 


0.69 juK 2 


2.38 juK 2 


L3 


8.4 


1.3 


0.9 


1.96 juK 2 


0.75 juK 2 


2.71 juK 2 



*L1 = 163Mpc_g8.7_130S = HMACH+LMACH, early reionization; 
L2 = 163Mpc_gl.7_8.7S = HMACH+LMACH, late reionization; 
L2M1 = 163Mpc_gl.7_8.7S_M300_J0.1 = HMACH+LMACH+MH; 
L3 = 163Mpc_g21.7_0 = HMACH only, early reionization. 



^Post-reionization kSZ values are CSF model of [59] for baryon correction. 



results for rms fluctuations and skewness of the maps smoothed with beamwidth and 
bandwidth like those for the LOFAR survey, for the two cases, which share the same 
global mean ionization history, finding good agreement. The fluctuation power spectrum 
is compared in Figure 1 (right panel) for the frequency which corresponds to the epoch 
when the global mean ionized fraction was 0.77 for both cases, showing that the larger 
box has more power at small wavenumber, as expected from the finite box size effect 
on the small box. At higher k, however, the agreement is better. These results will be 
presented in more detail in Iliev et al. (2012, in prep). An illustrative sky map of the 
21cm brightness fluctuations from the big box simulation at observer frequency 1 15.827 
MHz from z = 11.26 early in EOR, for a LOFAR-like beam (3 arcmin and average 
signal is zero) is shown in Figure 2 (left panel). 

3.3. Which Galaxies Reionized the Universe? 

Here we use a series of 163 Mpc box simulations to ask if fluctuations in the radiation 
backgrounds (21cm, CMB — polarization and kinetic Sunyaev-Zel'dovich, near-IR) 
from the EOR will allow us to determine which of the three types of galactic sources 
discussed above contribute most significantly to reionization. To keep track of the 
different cases used to address this question, the reader is referred to Table 1 and the 
meaning of the parameter labels for different cases. 

Can 21cm fluctuations tell us which galaxies reionized the universe? As shown in 
Figure 2 (right panel), for an illustrative epoch when global ionized fraction is x ( - = 
0.77, the difference between the curves for the case with HMACHs-only and those for 
the cases with LMACHS, with and without MHs, is larger than the error bars for a 
measurement of the power spectrum by a survey like MWA, so we might be able to 
distinguish the HMACH-only case, with higher power, from the cases with LMACHs, 
with or without MHs. But it cannot distinguish the case of HMACH + LMACH without 
MHs from that with MHs. 

Can CMB polarization fluctuations tell us which galaxies reionized the universe ? 
As shown in Figure 3 (left panel), of three cases that all finish ionizing at the same 
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FIGURE 2. 21cm brightness temperature fluctuations. (Left panel) Sky map for 607 Mpc box simula- 
tion observed at 1 15.827 MHz (i.e. looking back to z = 1 1 .26, early in EOR), smoothed for LOFAR-like 
beam (3' and average signal is zero). (Right panel) Comparison of power spectra at x, ,„ = 0.75 for 3 cases, 
as labelled. 



redshift ~ 8.3, the T es for the HMACH-only case is just below the WMAP7 1-sigma 
range, while both cases with LMACHs (with and without MHs) stay within that range, 
although the MH case has the highest T es . According to Figure 3 (right panel), even 
the higher sensitivity of Planck after 2-years of data will distinguish the HMACHs-only 
case from the others, but will hardly be able to distinguish the cases with HMACHs + 
LMACHs, with and without MHs, from each other, for such an early end to reionization. 
However, as Figure 4 (left panel) shows us, if reionization ended as late as z < 7, as some 
other observational evidence suggests, then the case with HMACHs + LMACHs (no 
MHs) makes T e s too small to be within the the WMAP7 1-sigma range, but HMACHs 
+ LMACHs + MHs is within this range. A more careful distinction requires the Planck 
2-year data soon to be released, however, but Figure 4 (right panel) clearly shows that, 
even if it is within the current WMAP uncertainties, the boost to T es and the polarization 
fluctuations caused by MHs should be readily detectable by Planck if reionization ended 
this late. Planck will thereby see the signature of the first stars at high redshift, currently 
undetectable by any other probe [48] . 

Can near-IR background fluctuations tell us which galaxies reionized the universe? 
The same stars that supply the ionizing radiation that escapes from galactic sources of 
reionization emit radiation below the H ionization threshold which escapes directly from 
the galaxies (in the absence of dust), and the absorbed fraction (1 — / esc ) of their ionizing 
photons is reprocessed into nebular emission that also escapes. This makes reionization 
sources also sources of the near-IR background, and their spatial clustering and its 
evolution contribute angular fluctuations in this background which probes the mass 
range of halos responsible for reionization, as well. We have applied our simulations 
of reionization to predict this fluctuating background in the near-IR, previously without 




FIGURE 3. (Left, top) Global ionization histories for three cases. (Left, bottom) Electron-scattering 
optical depths for three cases, compared with value from CMB polarization fluctuations; (Right) Predicted 
CMB polarization fluctuation angular power spectra for same three cases. Error bars are estimated Planck 
2-yr, 1-C7 sensitivity including cosmic variance. 




FIGURE 4. Same as Figure 3, except for different cases, including both early and late reionization 
cases. 



the LMACHs and their self -regulation [53], and now with them [54]. The talk by Beth 
Fernandez described this work in more detail. The angular power spectrum C\ tends to 
be dominated by the galaxies responsible for completing reionization (e.g. z ~ 6). The 
shape of the angular power spectrum C\ at high Z is sensitive to the amount of nonlinear 
bias of these galaxies relative to the total matter density, and since this bias depends on 
halo mass, the shape can be used to tell us what galactic halo mass was responsible for 
completing reionization. As Figure 5 shows, we find that Q is steeper for the case with 
suppression of LMACHs than with LMACHs included without suppression. In all cases, 
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FIGURE 5. (Left) Angular power spectra of the NIRB fluctuations for sources assumed to be Pop II 
galaxies with / e sc =0.1 and Salpeter IMF (after subtracting the shot-noise contribution), normalized to 
same amplitude at / w 5000; (Right) Same model, predictions in various bands compared to observational 
results. From [54], 
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FIGURE 6. Angular fluctuation power spectra for the kSZ effect, including contributions from post- 
reionization and the EOR for several reionization cases, as labelled (see Table 1), along with the SPT 
upper limits from [55]. From Park et al. 2012, in prep. 



we do not see a turn-over toward high I in the shape of I Q. 

Can the kinetic Sunyaev-Zel 'dovich effect from patchy reionization tell us which 
galaxies reionized the universe? The kSZ effect is a CMB temperature anistropy in- 
duced by electron scattering by free electrons moving along the line-of-sight, and is dis- 
tinguished from the related thermal SZ effect induced by the scattering by hot electrons 
in the intracluster gas, by the spectral difference between these two effects. As an inte- 
grated effect over the path of photons from recombination to the present, the kSZ signal 
is the sum of the contributions from the IGM during the EOR and post-reionization era. 
Recent results by the South Pole Telescope (SPT) have detected the tSZ effect and placed 
an upper limit on the kSZ effect at arcminute scale, I ~ 3000 [55], which has been inter- 
preted as a constraint on the duration and timing of reionization, imposed by the need to 
keep the EOR contribution below the upper limits, after subtracting the post-reionization 



signal expected theoretically, from the upper limit on the total kSZ signal [56, 57]. Non- 
detection by the SPT yields an angular fluctuation power upper limit on the total kSZ 
signal at / = 3000 of D/ < 2.8juK 2 (95% confidence), where D z = 1(1+1)Q/2jc, but if 
there is a correlation between galactic emission at these wavelengths and the tSZ effect, 
this introduces an uncertainty which raises the upper bound to the more conservative 
level of 6fiK 2 [55]. 

Our earlier predictions of the EOR contribution to this kSZ without LMACHs [40, 
58], have now been replaced by our new 163 Mpc box simulations of self-regulated 
reionization, listed in Table 1 . The results are summarized in Table 1 and Figure 6, from 
Park et al. (2012, in prep). All cases are easily consistent with the more conservative 
bound if there is a tSZ-CIB correlation, but they are also all allowed by the tighter 
bound, too. Apparently, the global reionization histories for these models include a 
range of durations which are larger than were found to be allowed by comparison of 
semi-numerical reionization models with the SPT upper limits by [56] and [57]. Alas, 
it does not yet seem possible to use the SPT limits on the kSZ fluctuations from patchy 
reionization to determine which galaxies reionized the universe, even though the results 
are beginning to place a significant constraint on the history of reionization. 

ACKNOWLEDGMENTS 

This study was supported in part by U.S. NSF grants AST-0708176 and AST- 1009799; 
NASA grants NNX07AH09G, NNG04G177G, and NNX1 1AE09G; and Chandra grant 
SAO TM8-9009X, and by the Swedish Research Council grant 2009-4088. KA. was 
supported in part by NRF grant funded by the Korean government MEST (Nos. 2009- 
0068141, 2009- 0076868, 2012R1A1A1014646, and 2012M4A2026720). I.T.I, was 
supported by The Southeast Physics Network (SEPNet) and the Science and Technology 
Facilities Council grants ST/F002858/1 and ST/1000976/1. The authors acknowledge the 
TeraGrid and the Texas Advanced Computing Center (TACC) at The University of Texas 
at Austin (URL: http://www.tacc.utexas.edu), and the Swedish National Infrastructure 
for Computing (SNIC) resources at HPC2N (UmeN, Sweden) for providing HPC and 
visualization resources that have contributed to the results reported within this paper. 

REFERENCES 

1. P. R. Shapiro, and M. L. Giroux, ApJ 321, L107-L1 12 (1987). 

2. P. R. Shapiro, M. L. Giroux, and A. Babul, ApJ 427, 25-50 (1994). 

3. N. Y. Gnedin, ApJ 542, 535-541 (2000), arXiv:astro-ph/00021 51 . 

4. P. R. Shapiro, I. T. Iliev, and A. C. Raga, MNRAS 348, 753-782 (2004), arXiv:astro-ph/0307266. 

5. Z. Haiman, T. Abel, and M. J. Rees, ApJ 534, 1 1-24 (2000), arXiv:astro-ph/9903336. 

6. H. Susa, PASJ 58, 445^60 (2006), arXiv:astro-ph/0601 642. 

7. J. H. Wise, and T. Abel, MNRAS 414, 3458-3491 (201 1), 1 01 2.2865. 

8. A. H. Pawlik, and J. Schaye, MNRAS 389, 651-677 (2008), 0802.1 71 5. 

9. M. Petkova, and V. Springel, MNRAS 412, 935-946 (201 1), 1 008.4459. 

10. B. Semelin, F. Combes, and S. Baek, A&A 474, 365-374 (2007), 0707.2483. 

11. T. Okamoto, K. Yoshikawa, and M. Umemura, MNRAS 419, 2855-2866 (2012), 1 1 05.3269. 

12. S. Cantalupo, and C. Porciani, MNRAS 411, 1678-1694 (201 1), 1 009.1 625. 

13. M. Petkova, and U. Maio, MNRAS 422, 3067-3080 (2012), 1 1 1 0.0492. 



14. M. Petkova, and V. Springel, MNRAS 396, 1383-1403 (2009), 081 2.1 801 . 

15. K. Finlator, F. Ozel, and R. Dave, MNRAS 393, 1090-1 106 (2009), 0808.3578. 

16. G. Altay, R. A. C. Croft, and I. Pelupessy, MNRAS 386, 1931-1946 (2008), 0802.3698. 

17. R. A. C. Croft, and G. Altay, MNRAS 388, 1501-1520 (2008), 0709.2362. 

18. H. Y. Trac, and N. Y. Gnedin, Advanced Science Letters 4, 228-243 (201 1), 0906.4348. 

19. D. Aubert, and R. Teyssier, ApJ 724, 244-266 (2010), 1 004.2503. 

20. D. Aubert, and R. Teyssier, MNRAS 387, 295-307 (2008), 0709. 1 544. 

21. S. Baek, P. Di Matteo, B. Semelin, F. Combes, and Y. Revaz, A&A 495, 389-405 (2009), 0808.0925. 

22. H. Trac, and R. Cen, ApJ 671, 1-13 (2007), arXiv:astro-ph/061 2406. 

23. M.-S. Shin, H. Trac, and R. Cen, ApJ 681, 756-770 (2008), 0708.2425. 

24. M. McQuinn, A. Lidz, O. Zahn, S. Dutta, L. Hernquist, and M. Zaldarriaga, MNRAS 377, 1043-1063 
(2007), arXiv:astro-ph/061 0094. 

25. M. A. Alvarez, V. Bromm, and P. R. Shapiro, ApJ 639, 621-632 (2006), arXiv:astro-ph/0507684. 

26. K. Kohler, N. Y Gnedin, J. Miralda-Escude, and P. A. Shaver, ApJ 633, 552-559 (2005), arXiv: 
astro-ph/0501086. 

27. A. Maselli, B. Ciardi, and A. Kanekar, MNRAS 393, 171-178 (2009), 081 1 .031 5. 

28. A. Sokasian, T Abel, L. Hernquist, and V. Springel, MNRAS 344, 607-624 (2003), arXiv:astro-ph/ 
0303098. 

29. A. O. Razoumov, M. L. Norman, T Abel, and D. Scott, ApJ 572, 695-704 (2002), arXiv:astro-ph/ 
0109111. 

30. A. Maselli, A. Ferrara, and B. Ciardi, MNRAS 345, 379-394 (2003), arXiv:astro-ph/03071 1 7. 

31. N. Y. Gnedin, and T. Abel, New Astronomy 6, 437-455 (2001), arXiv:astro-ph/01 06278. 

32. T. Abel, M. L. Norman, and P. Madau, ApJ 523, 66-71 (1999), arXiv:astro-ph/981 2151. 

33. B. Ciardi, J. S. Bolton, A. Maselli, and L. Graziani, MNRAS 423, 558-574 (2012). 

34. B. Ciardi, F. Stoehr, and S. D. M. White, MNRAS 343, 1 101-1 109 (2003), arXiv:astro-ph/0301 293. 

35. N. Y. Gnedin, ApJ 535, 530-554 (2000), arXiv:astro-ph/9909383. 

36. I. T. Iliev, B. Ciardi, M. A. Alvarez, A. Maselli, A. Ferrara, N. Y. Gnedin, G. Mellema, T. Nakamoto, 
M. L. Norman, A. O. Razoumov, E.-J. Rijkhorst, J. Ritzerveld, P. R. Shapiro, H. Susa, M. Umemura, 
and D. J. Whalen, MNRAS 371, 1057-1086 (2006), arXiv:astro-ph/06031 99. 

37. I. T. Iliev, D. Whalen, G. Mellema, K. Ahn, S. Baek, N. Y. Gnedin, A. V. Kravtsov, M. Norman, 
M. Raicevic, D. R. Reynolds, D. Sato, P. R. Shapiro, B. Semelin, J. Smidt, H. Susa, T. Theuns, and 
M. Umemura, MNRAS 400, 1283-1316 (2009), 0905.2920. 

38. I. T. Iliev, G. Mellema, U.-L. Pen, H. Merz, P. R. Shapiro, and M. A. Alvarez, MNRAS 369, 1625- 
1638 (2006), arXiv:astro-ph/051 21 87. 

39. H. Merz, U.-L. Pen, and H. Trac, New Astronomy 10, 393-407 (2005), arXiv:astro-ph/0402443. 

40. I. T. Iliev, G. Mellema, U.-L. Pen, J. R. Bond, and P. R. Shapiro, MNRAS 384, 863-874 (2008), 
arXiv:astro-ph/0702099. 

41. I. T. Iliev, G. Mellema, P. R. Shapiro, and U.-L. Pen, MNRAS 376, 534-548 (2007), arXiv:astro-ph/ 
0607517. 

42. I. T. Iliev, P. R. Shapiro, G. Mellema, H. Merz, and U.-L. Pen, "Simulating Cosmic 
Reionization," in "Proceedings of the TeraGrid 2008 Conference. June 9-13, 2008. 
Las Vegas, USA. <A href="http://archive.teragrid.org/events/teragrid08/index.htm"> 
http://archive.teragrid.0rg/events/teragridO8/index.htnK/A>, p. 3 1 " ', 2008, p. 31, 0806.2887. 

43. P. R. Shapiro, I. T. Iliev, G. Mellema, U.-L. Pen, and H. Merz, "The Theory and Simulation of the 21- 
cm Background from the Epoch of Reionization," in The Evolution of Galaxies Through the Neutral 
Hydrogen Window, edited by R. Minchin, and E. Momjian, 2008, vol. 1035 of American Institute of 
Physics Conference Series, pp. 68-74, 0806.3091 . 

44. I. T. Iliev, G. Mellema, P. R. Shapiro, U.-L. Pen, Y. Mao, J. Koda, and K. Ahn, MNRAS 423, 2222- 
2253(2012), 1107.4772. 

45. M. M. Friedrich, G. Mellema, I. T. Iliev, and P. R. Shapiro, MNRAS 421, 2232-2250 (2012), 
1201.0602. 

46. K. K. Datta, M. M. Friedrich, G. Mellema, I. T. Iliev, and P. R. Shapiro, MNRAS 424, 762-778 
(2012), 1203.0517. 

47. K. Ahn, P. R. Shapiro, I. T. Iliev, G. Mellema, and U.-L. Pen, ApJ 695, 1430-1445 (2009), 0807. 
2254. 

48. K. Ahn, I. T. Iliev, P. R. Shapiro, G. Mellema, J. Koda, and Y. Mao, ApJ 756, L16 (2012), 1 206.5007. 



49. H. Jensen, P. Laursen, G. Mellema, I. T. Iliev, J. Sommer-Larsen, and P. R. Shapiro, ArXiv e-prints 
(2012), 1206.4028. 

50. Y. Mao, P. R. Shapiro, G. Mellema, I. T. Iliev, J. Koda, and K. Ahn, MNRAS 422, 926-954 (2012), 
1104.2094. 

51. K. K. Datta, G. Mellema, Y. Mao, I. T. Iliev, P. R. Shapiro, and K. Ahn, MNRAS 424, 1877-1891 
(2012), 1109.1284. 

52. M. M. Friedrich, G. Mellema, M. A. Alvarez, P. R. Shapiro, and I. T. Iliev, MNRAS 413, 1353-1372 
(2011), 1006.2016. 

53. E. R. Fernandez, E. Komatsu, I. T. Iliev, and P. R. Shapiro, ApJ 710, 1089-1 1 10 (2010), 0906.4552. 

54. E. R. Fernandez, I. T. Iliev, E. Komatsu, and P. R. Shapiro, ApJ 750, 20 (2012), 1 1 1 2.2064. 

55. C. L. Reichardt, L. Shaw, O. Zahn, K. A. Aird, B. A. Benson, L. E. Bleem, J. E. Carlstrom, C. L. 
Chang, H. M. Cho, T. M. Crawford, A. T. Crites, T. de Haan, M. A. Dobbs, J. Dudley, E. M. 
George, N. W. Halverson, G. P. Holder, W. L. Holzapfel, S. Hoover, Z. Hou, J. D. Hrubes, M. Joy, 
R. Keisler, L. Knox, A. T. Lee, E. M. Leitch, M. Lueker, D. Luong-Van, J. J. McMahon, J. Mehl, 
S. S. Meyer, M. Millea, J. J. Mohr, T. E. Montroy, T. Natoli, S. Padin, T. Plagge, C. Pryke, J. E. Ruhl, 
K. K. Schaffer, E. Shirokoff, H. G. Spieler, Z. Staniszewski, A. A. Stark, K. Story, A. van Engelen, 
K. Vanderlinde, J. D. Vieira, and R. Williamson, ApJ 755, 70 (2012), 1111 .0932. 

56. O. Zahn, C. L. Reichardt, L. Shaw, A. Lidz, K. A. Aird, B. A. Benson, L. E. Bleem, J. E. Carlstrom, 
C. L. Chang, H. M. Cho, T. M. Crawford, A. T. Crites, T. de Haan, M. A. Dobbs, O. Dore, J. Dudley, 
E. M. George, N. W. Halverson, G. P. Holder, W. L. Holzapfel, S. Hoover, Z. Hou, J. D. Hrubes, 
M. Joy, R. Keisler, L. Knox, A. T. Lee, E. M. Leitch, M. Lueker, D. Luong-Van, J. J. McMahon, 
J. Mehl, S. S. Meyer, M. Millea, J. J. Mohr, T. E. Montroy, T. Natoli, S. Padin, T. Plagge, C. Pryke, 
J. E. Ruhl, K. K. Schaffer, E. Shirokoff, H. G. Spieler, Z. Staniszewski, A. A. Stark, K. Story, A. van 
Engelen, K. Vanderlinde, J. D. Vieira, and R. Williamson, ApJ 756, 65 (2012), 1111 .6386. 

57. A. Mesinger, M. McQuinn, and D. N. Spergel, MNRAS 422, 1403-1417 (2012), 1 1 1 2.1 820. 

58. I. T. Iliev, U.-L. Pen, J. R. Bond, G. Mellema, and P. R. Shapiro, ApJ 660, 933-944 (2007), 
arXiv:astro-ph/0609592. 

59. L. D. Shaw, D. H. Rudd, and D. Nagai, ApJ 756, 15 (2012), 1 1 09.0553. 



